Non-equilibrium dynamics of the driven Hubbard model. 
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We investigate the dynamics of a two-dimensional Hubbard model in a static electric field in order 
to identify the conditions to reach a non-equilibrium stationary state. For a generic electric field, 
the convergence to a stationary state requires the coupling to a thermostating bath absorbing the 
work done by the external force. Following the real-time dynamics of the system, we show that 
a non-equilibrium stationary state is reached for essentially any value of the coupling to the bath. 
We map out a phase diagram in terms of dissipation and electric field strengths and identify the 
dissipation values in which steady current is largest for a given field. 
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The theoretical investigation of out-of-equilibrium 
strongly correlated quantum systems recently generated 
tremendous interest, stimulated by novel experimen- 
tal techniques which allow to explore transport proper- 
ties of correlated materials in a non-equilibrium regime, 
notwithstanding the recent achievements in the field of 
optically trapped cold-atoms. These experimental ad- 
vances posed serious challenges to the theory, which re- 
quire the development of new ideas to be tackled. In 
fact, although we achieved a fair understanding of corre- 
lated materials at equilibrium, we are just beginning to 
uncover the basic principles governing quantum systems 
far from equilibrium. 

Dynamical mean-field theory (DMFT) is an estab- 
lished method to investigate correlated materials [1]. The 
recent extension of the DMFT out of equilibrium [2, 3], 
provides us with a reliable tool to clarify how correla- 
tion effects influence the non-equilibrium dynamics of 
quantum systems. Non-equilibrium DMFT has been suc- 
cessfully applied to study quantum quenches - sudden 
changes of some control parameter [4] , and to investigate 
driven correlated systems [2, 5]. In this context, DMFT 
has been used to show how interactions favor the forma- 
tion of stationary states by suppressing Bloch oscillations 
of the currcnt[2. 5, 6] and to analyze the dielectric break- 
down of Mott insulators [7]. 

In this work we focus on the role of the coupling to an 
external thermostat in the out-of-equilibrium dynamics 
of a correlated system and we show when and how it leads 
to a non-equilibrium stationary state (NSS). Despite its 
importance, the role of dissipation in driven correlated 
systems has only been discussed assuming the existence 
of a NSS [8], while we are not aware of any study following 
the real-time dynamics leading to such state. In the clas- 
sical framework, dissipation is usually introduced by cou- 
pling the system to a set of reservoirs that ultimately im- 
pose a suitable constraint on the equations of motion [9]. 
Nevertheless, the direct extension of this approach to the 



quantum regime is not straightforward because of the 
Hamiltonian nature of the quantum equations[10]. 

To include dissipation effects we couple the Hubbard 
model, the paradigm of strong correlations, to a fermionic 
bath. We follow the real-time dynamics when the system 
is driven out of equilibrium by a static electric field. We 
show that, for a given value of the electric field, coupling 
to the thermostat is necessary to reach a physically rele- 
vant NSS. Remarkably, a NSS with finite current can be 
reached for almost any value of the dissipation regardless 
the initial conditions. The coupling to a bath is therefore 
also an essentially sufficient condition to dynamically ap- 
proach the NSS. Finally, we characterize the properties 
of the NSS in terms of local observables and construct 
the phase-diagram of the model. 

Model. We consider the two-dimensional Hubbard 
model on a square lattice with spacing a. The model 
Hamiltonian reads {c ~ 1): 

'Hc = Y^ e{k - eA{t))cl^cua + ^i,T"*4 (1) 



where fiia = cJ^Cia — 1/2. This model describes a half- 
filled band of conduction electrons with dispersion £(k) 
and hopping amplitude — I (setting the energy unit 
of the problem) subject to a local Coulomb repulsion of 
strength U. The system is coupled to a constant and ho- 
mogeneous electric field E, derived from a purely vector 
potential A{t) = -r{t)m, with E and Q ^ (7r,7r). 

The function r{t) specifies the switching protocol of the 
field. In this work we consider r{t) = 6{t — to), cor- 
responding to a sudden switch on of the electric field. 
Different choices for r{t) will be explicitly stated. The 
electric field unit is eEa. The coupling to electric field is 
realized via the Peierls substitution: k — > k = k — eA(t). 

In a solid state system the work done by the electric 
field on the electrons is constantly transformed into heat 
by various scattering mechanisms. Thus, coupling to an 
external thermostat is crucial to maintain the internal en- 
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ergy conserved on average. A thermostating mechanism 
is realized by assuming the system to be locally coupled 
to a bath of non-interacting electrons at a temperature 
T ~ 0.01 and which are unaffected by the electric field. 
As the details of the internal structure of the thermostat- 
ing reservoirs are irrelevant with respect to the physics 
of the NSS, we consider a set of identical systems with 
constant density of states with a bandwidth W . The 
coupling has the form: 



For simplicity we choose Vk^i = V and wc drop any refer- 
ence to spin index. 

The Kadanoff-Baym-Kcldysh[ll] formalism provides a 
natural framework to investigate the dynamics in out-of- 
equilibrium conditions. In the presence of sizeable elec- 
tric fields the non-equilibrium dynamics is expected to 
be driven by the field and by the coupling to the thermo- 
stat. Being interested in the relaxation towards the NSS, 
we can simplify the treatment considering only the real- 
time branches in the Keldysh contour and neglecting the 
imaginary-time segment. As a consequence, at time to 
we will simultaneously switch-on the electric field E and 
the interaction J7, introducing also a correlation quench. 
This is not expected to strongly influence the results in 
the presence of the thermostat, which is able to dissipate 
the extra energy involved in the quench as well as the 
energy pumped by the field [12]. 

The thermostated system obeys the following Dyson 
equation on the Keldysh contour C 



Gk — Gok + Gqk ■ S„ • Gk 



(2) 



where all quantities represent continuous operator of two 
time variables {t, t') £ C, the symbol • denotes the convo- 
lution product f ■ g = dzf{t, z)g(z, t') and denotes 
the Keldysh self-energy function. The non-interacting 
lattice Green's function Qq,^ is given by: 



G^^{t,t')^[^dt-e{n)]-6c{tX) 



(3) 



where 5c indicates the delta function on the contour 
C, So(t - t') = V^git - t') and g{t - t') the Fourier 
transform of the non-interacting local bath Green's func- 
tion g{uj) = \[n\{W + uj)/{W ~u)\-iTTe{W -\uj\)]/2W. 
Thus, the effective coupling to the thermal reservoirs is 
A = V^/2W. 

Applying Q^^^ on both sides of the Eq. 2 and using 
Eq. 3, we recast the Dyson equation into the Kadanoff- 
Baym equation, that is the equation of motion for the 
Keldysh Green's function, which reads[13]: 



G<- 



e{K) 
-idt' 



G<-j:t 



(4) 
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Figure 1. (Color online) (a) dynamics of the local current J 
for U — 6, E — 1.26 and increasing coupling to thermostat 
A. (b) same as (a) and U = 4 and E = 4.7. (c) dynamics 
of the local current from different realizations of the initial 
conditions: electric field quench on the non-interacting U = 
(black) and interacting U — 4 state (red), smooth switch-on 
of the electric field for (7 = and r = 2 (blue) (see text). The 
inset show the profile of the electric field, (d) linear/non- 
linear crossover of the steady current as a function of the 
electric field E. Data for A = 0.025. 



As an effect of the symmetries relating the Keldysh 
components, Eqs. 4 determine a system of coupled ordi- 
nary first-order differential equations in the (t,t')-plane 
for the lattice Keldysh Green's function G^. The nu- 
merical solution of this systcm[14] permits to obtain Gk 
up to arbitrary large times, once suitable initial con- 
ditions arc specified. We fixed the time-step parame- 
ter, introduced by discretization of the contour C, to 
6t < 0.1. The initial conditions for the system 4 are 
given specifying the value of momentum distribution at 
an (arbitrary) initial time t = t^. In this Letter we set 
nk(to 0) = -iG<(0,0) = /(e(k)), with / the Fermi- 
Dirac distribution, except where otherwise stated. 

We use DMFT to deal with correlations in a non- 
perturbative way. Within DMFT the lattice self-energy 
is replaced by its local component, so that momentum 
dependence is determined by the non-interacting dis- 
persion. The self-energy can be derived from the self- 
consistent solution of an impurity problem written in 
terms of a local Weiss field Qo- This latter describes 
the effective medium coupled to the impurity and has 
to be self-consistently determined solving the following 
equations: 



t/o — i • ^\oc ■ 

— F • G^ • Ft 



F = [I- 



G\oc 



(5) 



where local functions arc obtained summing the solution 
of Eqs. 4 for any k G [~7r,7r]^: Fioc = Zlk-^«- To close 
the DMFT equations it is necessary to determine the 
impurity self-energy I](i, t'). In this work we use iterated 
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Figure 2. (Color online) Left panel: effective temperature 
evolution T^si^) for U = 6, E — 1.9. Right panel: effective 
temperature TcS as a function of the electric field, for U = 4 
and increasing coupling to thermostat. 




second-order perturbation theory in U[l]: 

j:^{t,n = u'[g^{t,t'rgf{t',t) (6) 

Equations 4, 5 and 6 are iterated until a self-consistent 
solution is obtained. 

Results. The approach to the stationary non- 
equilibrium state can be characterized following the real- 
time dynamics of suitable observables, such as the local 
current J(t) = —ie/7r VKG'^(t, t), where v„ = V„£(k) 
is the electronic velocity field. We focus on the correlated 
metallic phase. Our results for the local current are pre- 
sented in Fig. 1. The application of a constant electric 
field on a periodic lattice structure produces an oscillat- 
ing current (Bloch oscillations) of frequency lob = eEa. 
In the absence of coupling to the external bath (A = 0), 
the electron-electron interactions suppress the Bloch os- 
cillations leading to an exponentially decaying current 
which eventually converges to zero at very long time [see 
panels (a)-(b)]. However, the relaxation process shows 
two different regimes depending on the value of the in- 
teraction U as found in [6]. 

The non-equilibrium dynamics of the system dramati- 
cally changes in presence of finite coupling to the thermo- 
stat (A > 0). The local current relaxes to a finite value 
corresponding to the formation of a NSS. This effect is de- 
tailed in panels (a)-(b) for two different values of the elec- 
tric field. The approach to the NSS is found to be inde- 
pendent from the initial conditions, confirming that the 
non-equilibrium physics is governed by the field and the 
dissipation term. In panel (c) of Fig. 1 we compare the 
solutions obtained using initial conditions r{t) ~ 6{t) and 
'^k(O) = /(£(k)) with the solutions obtained: (i) starting 
from the interacting equilibrium momentum distribution 
for [/ = 4, (ii) using smooth switching of the electric 
field r{t) = [1 - 3/2cos(7ri/T) + l/2cos{T:t/Tf]/2. The 
convergence to the same NSS is evident. 

In panel (d) of Fig. 1 we plot the steady current which 
shows a linear/nonlinear crossover as a function of the 
electric field [cf. panel (d)][15, 16]. At small fields the 
current is linear in E, as expected by continuity with 
perturbed equilibrium state, then it reaches a maximum 
before decreasing as the field is further increased. 



Figure 3. (Color online) (a) Non-equilibrium dynamics of the 
double occupations (d) for U = 6, E = 1.25 and increasing A. 
(b) Double occupancy of the NSS as a function of the electric 
field E for U — 6. The equilibrium solution for the same value 
of the interaction (dotted line) is reported for comparison, (c) 
Entropy 5* (see text) as a function of the electric field E, for 
[7 = 4 and the same value of A as in (b). 

A better understanding of the effects introduced by 
inclusion of a thermostat is obtained looking at the effec- 
tive temperature Teg. This is defined as the temperature 
associated to an equilibrium solution with the same en- 
ergy V,{t) = (K) + (V) of the non-equilibrium state[6, 7] 
and the same value of the interaction. For A = the ef- 
fective temperature rapidly diverges due to the Coulomb 
scattering between accelerated electrons. In the presence 
of dissipation the effective temperature relaxes instead 
to a finite value. The coupling to the thermostat pro- 
vides a competing scattering channel which prevents the 
population of high-momentum states and allows the con- 
servation (on average) of the total energy and a constant 
effective temperature. 

In the right panel of Fig. 2 we show the behavior of 
the effective temperature of the NSS. For any given cou- 
pling A, Teff is a monotonically increasing function of the 
electric field, so that arbitrarily large fields imply arbi- 
trarily large T^g. For large fields T^s is depressed by the 
increasing of the coupling A. [However at smaller fields 

< 1 we observe a non-monotonic behavior of the effec- 
tive temperature, due the increased scattering with bath 
degrees of freedom.] 

A key quantity to characterize correlated systems is 
the double occupancy (d) which is proportional to the po- 
tential energy and measures the effectiveness of correla- 
tions. In panel (a) of Fig. 2 we illustrates the formation of 
NSS from the dynamics of (d). In the non-thermostated 
case (A = 0) (d) relaxes to the non-interacting value 
(d) = 1/4. Conversely, a finite coupling to the ther- 
mostat induces the relaxation to a NSS, which is char- 
acterized by a smaller value (d) < 1/4. The stationary 
double occupancy, plotted in panel (b) as a function of 
E, is usually different from the corresponding equilibrium 
value. Remarkably, the evolution of (d) as a function of 
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Figure 4. (Color online) Top panel: phase-diagram of the 
thermostated non-equilibrium Hubbard model. Diagram is 
obtained from the normalized local current J as a function 
of electric field strength E and coupling to thermal bath A. 
Horizontal lines indicates the cuts shown in the bottom panel. 
Bottom panel: local current J and entropy S behavior as a 
function of A, for E = 1 (circles), E = 2 (triangles), E = 3 
(squares). 



the electric field presents a minimum at a characteris- 
tic value of the field strength. This behavior suggests 
that small applied fields, in conjunction with finite dis- 
sipation, increase the degree of electron localization in 
the NSS. Since the effective temperature of the NSS in- 
creases with the electric field, we can connect this result 
with the Pomeranchuk effect observed in equilibrium. In 
the latter case the more localized state is favored at finite 
temperature because of its larger spin entropy [1]. Upon 
increasing A, {d) rapidly reduces its range of variation 
getting closer to the equilibrium value (A = E = 0) for 
all values of i?, implying that, in terms of correlation 
properties, the NSS becomes closer to the equilibrium 
state as the dissipation increases. 

Further information about the NSS can also be ob- 
tained evaluating the function S ~ — 2^^nKln(nK) 
which would coincide with the equilibrium entropy at 
U = 0. The extension of entropy for non-equilibrium sys- 
tems is a debated issuc[17] which goes beyond the aim 
of this paper. Therefore we use S simply as a tool to 
extract information about the approach to a stationary 
state and we do not interpret it as a physical entropy. In 
the absence of a thermostat 5* dynamically saturates to 
its high-temperature limit S — In 4, corresponding to the 
approach to the zero current stationary state. A substan- 
tial reduction of S in the non-linear regime is achieved 
by coupling the system to the thermostat, as detailed in 
Fig. 2(c). Nevertheless, for a given coupling A, S in the 
NSS is a monotonically increasing function of the electric 
field and reaches its limiting value for large enough elec- 
tric fields, in agreement with the unbounded increase of 
the effective temperature at constant coupling A. 

The dependence of the steady current on the coupling 



A and the electric field can be cast in a diagram in the 
{E-A) plane shown in Fig. 4. At the borders of the di- 
agram we identify two small regions where the systems 
shows small or null steady current. Near the A = axis, 
the energy injected by the field largely overcomes that 
absorbed by the thermostat, so that the continuous heat- 
ing ultimately leads to a reduction of the steady current. 
Conversely, near the E = axis, the large scattering in- 
troduced by strong coupling to thermostat reduces the 
linear conductivity and thus the corresponding steady 
current. In a wide intermediate region centered around 
the diagonal of the phase-diagram {E/2tt ~ A) we found 
largest values of the steady current. In this region the dis- 
sipation is sufficient to get rid of the extra energy pumped 
in by the field, but it is not too large to overcome the ef- 
fect of the field driving a coherent current. The bottom 
panel shows the current and S along the three horizontal 
cuts in the diagram, underlining that in the intermediate 
region the current is maximal and S is minimal. 

Conclusions. Using DMFT in combination with a di- 
rect solution of the Kadanoff-Baym equations we in- 
vestigated the non-equilibrium dynamics of the two- 
dimensional driven Hubbard model coupled to electronic 
reservoirs. We reported that for a generic value of the 
field the coupling to an external bath is a necessary and, 
remarkably, also a sufficient condition to reach a non- 
stationary steady state with a finite current. We charac- 
terized the properties of the NSS in terms of experimen- 
tally accessible quantities and studied their dependence 
on the coupling to the thermostat, identifying the con- 
ditions to obtain a maximum of the steady current for a 
given field. Our work is a fundamental step towards a sat- 
isfactory description of non-equilibrium solids, in which 
a certain degree of dissipation is always present. An ex- 
plicit coupling to a thermostat is shown to be essential 
to obtain a description of a non-equilibrium stationary 
states with finite effective temperature, in contrast to a 
modelling which neglects dissipation effects and that can 
only give rise to transient states whose relevance for the 
physics of actual materials remains questionable. 
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